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Abstract. 

We propose a physical model for the nonlinear inelastic mechanics of sticky 
biopolymer networks with potential applications to inelastic cell mechanics. It 
consists in a minimal extension of the glassy wormlike chain (GwLC) model, which 
has recently been highly successful as a quantitative mathematical description of 
the viscoelastic properties of biopolymer networks and cells. To extend its scope 
to nonequilibrium situations, where the thermodynamic state variables may evolve 
dynamically, the GwLC is furnished with an explicit representation of the kinetics 
of breaking and reforming sticky bonds. In spite of its simplicity the model exhibits 
many experimentally established non-trivial features such as power-law rheology, stress 
stiffening, fluidization, and cyclic softening effects. 
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1. Introduction 

In many studies of cell mechanics and dynamics, the cell is characterized as a viscoelastic 
body [IHS]. It is an intriguing question to what extent such mechanical behaviour 
can be rationalized in terms of comparatively simple polymer physics models. In this 
respect, the comparison of cell rheological data and minimalistic in vitro reconstituted 
constructs of the cytoskeleton, such as pure actin solutions |1H8] or crosslinked actin 
networks [OHT^. has recently provided many new insights. Important progress has 
also been achieved in the development of phenomenological mathematical descriptions. 
This includes approaches related to the tube model [T5HT8], tensegrity-based approaches 
[T2H2T], effective-medium models [22H21], and some others [221 12S]- In particular, the 
glassy wormlike chain (GwLc) model [27], a phenomenological extension of the standard 
Wlc model of semifiexible polymers [2H] has been successful in describing a plethora 
of rheological data for polymer solutions |71 129] and living cells [30] over many decades 
in time with a minimum of parameters. However, all these studies were primarily 
concerned with viscoelastic behaviour, while the latest investigations have underscored 
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the glassy |3TH33] fragile [Ml |35], and inelastic [MH36] character of the mechanical 
response of living cells. Even for biopolymer networks in vitro, experiments operating 
in the nonlinear regime had so far to resort to special protocols that minimize plastic 
flow [H [37] in order to make contact with dedicated theoretical models. 

The aim of the present contribution is to overcome this restriction by extending 
the GwLC to situations involving inelastic deformations. As a first step, we concentrate 
onto reversible inelastic behaviour, where the deformation does not alter the microscopic 
ground state. The protocol applied by Trepat et al. [34j provides a paradigmatic 
example. Cells are subjected to a transient stretch such that, after some additional 
waiting time in the unstretched state, the (linear) material properties of the initial 
state are recovered. The simplification for the theoretical modelling results from the 
assumption that not only the macro-state but also the micro-state of the system may 
to a good approximation be treated as reversible under such conditions; i.e., we assume 
that the complete conformation of the polymer network, including the transiently 
broken bonds between adjacent polymers, is constrained to eventually return to its 
original equilibrium state. For the time-delayed hysteretic response of the network to 
such protocols one could thus still speak of a viscoelastic ( "anelastic" ) response in an 
operational sense, but we refrain from doing so in view of the fundamentally inelastic 
nature of the underlying stochastic process — in contrast to the reversible softening 
effects observed in [38], for example. Indeed, by simply allowing bonds to reform 
in new conformational states, the model developed below can readily be extended 
to arbitrary irreversible plastic deformations, as will be demonstrated elsewhere [15] . 
Before entering the discussion of our model, we would also like to point out that the 
proposed (inelastic) extension of the GwLC is strongly constrained by identifying the 
newly introduced parameters with those of the original (viscoelastic) GwLC model, 
where possible. Despite its increased complexity, the extended model will therefore 
enable us to subject the underlying physical picture to a more stringent test than 
hitherto possible by comparing its predictions to dedicated experiments. Moreover, 
unlike current state-of-the-art simulation studies [2S] it is not limited to rod networks 
but is firmly routed in a faithful mathematical description of the underlying Brownian 
polymer dynamics. 

This paper is organized as follows. First, we review some basic facts about the 



GwLC in section 2.1 Next, in section 2.2, we introduce our extended reversible inelastic 



version, which we formulate using the notion of an effective interaction potential as in the 
original construction of the G WLC in [27] . (A preliminary account of the basic procedure 
and some of its cell-biological motivation including reversible bond-breaking kinetics has 
recently been given in a conference proceedings [39].) Sections 3.1 and 3.2 explain the 
physical mechanism underlying the mechanical response under pulsed and periodically 
pulsed loading, while section 3^ illustrates its phenomenology. We demonstrate that 
the model exhibits the hallmarks of nonlinear cell mechanics: strain/stress stiffening, 
fiuidization, and cyclic softening [51H5B] . Section 3.4 investigates the relevance of the 
lately quantified structural heterogeneities in networks of semifiexible polymers [40j for 
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the mechanical properties, before we conclude and close with a brief outlook. 
2. Theory 

2.1. The glassy wormlike chain 

The glassy wormlike chain (GwLc) is a phenomenological extension of the wormlike 
chain (Wlc) model, the well-established standard model of semiflexible polymers. A 
broad overview over Wlc and GwLC dynamics can be found elsewhere [28]. The 
Wlc describes the mechanics of an isolated semiflexible polymer in an isothermal 
viscous solvent. In the weakly bending rod approximation, a solution of the stochastic 
differential equations of motion for the Wlc is possible via a mode decomposition ansatz 
for the transverse displacement of the polymer contour from the straight ground state. 
The individual modes labelled by an index n are independent of each other and decay 
exponentially with rates r~^. For convenience, we set the thermal energy ksT = 1, 
so that the bending rigidity can be identified with the persistence length ip, in the 
following. Using this convention, the Wlc expression for the transverse susceptibility 
of a polymer of contour length L (with respect to a point force) reads [27] 



Here, / is an optional backbone tension, and = ir'^ip/ the Euler buckling force for 
a hinged rod of length L. The different powers of n in the denominator give notice of 
the competition of bending and stretching forces. 



Figure 1. Cartoon of a test polymer (differently coloured) in a network. The trapping 
of the test polymer is schematically represented by means of an effective interaction 
(figure [2| indicated by the potential wells at sticky entanglement points, which are 
on average separated by an contour length Lg. The test polymer can bind/unbind by 
overcoming an energy barrier of height £. The average contour length between the 
closed bonds is A. 



In the GwLC, interactions of the polymer with the surrounding network (e.g. 
excluded volume interactions or sticky contacts) are reflected in an altered Wlc mode 




(1) 
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relaxation spectrum |27l |29]. The intuitive — albeit not fully microscopic — picture 
underlying the formulation of the GwLC is illustrated in figure [T| which depicts a 
test polymer reversibly bound to the background network via potential wells at all 
topological contacts, the so-called entanglement points. Consider now a generic point 
somewhere along the polymer backbone (not an entanglement point). It can relax 
freely until the constrictions are being felt, which slow down the contributions from 
long wavelength bending modes. The GwLC translates this intuition into a simple 
prescription for the mode spectrum: the short-wavelength modes are directly taken over 
from the wormlike chain model, while modes of wavelength A„ longer than the typical 
contour length A between adjacent bonds are modified to account for the slowdown. 
Motivated by the physical picture illustrated in figure [1} the slowdown of the relaxation 
of a wavelength A„ is expressed by an Arrhenius factor exp (A„/A — 1)] for breaking 
(An/A — 1) potential energy barriers of height S simultaneously. Accordingly, the 
phenomenological recipe to turn a Wlc into a GwLC reads: 

r ^ f = I ^"^ (2) 

" " \ r„exp[^(A„,/A-l)] A„ > A ' ^ ^ 

Upon inserting this into ([T| one obtains the GwLC susceptibility q;gwlc(i^)- The 
microscopic "modulus" for transverse point excitations of a generic backbone element 
on a test polymer is then defined as the inverse susceptibility, gcwLci^) = l/aGwLc(i^)- 
An approximate expression for the macroscopic shear modulus is obtained along similar 
lines |22l|27|. 

In the original equilibrium GwLC theory, A was assumed to be a constant on the 
order of the entanglement length of the network, A > Lg. Note, however, that Le 
is a geometric quantity (which is determined by the polymer concentration and the 
persistence length) while the contour length A between closed bonds clearly depends on 
the state of the bond network. One therefore has to allow for an increase of A with the 
number of broken bonds in non-equilibrium applications. This issue is explored in the 
following. 



2.2. Effective interaction potential and bond kinetics 

All mechanical quantities calculated within the GwLC model crucially depend on the 
interaction length A. In previous applications of the model [7l|271[30] it was assumed that 
A remains constant — equal to its equilibrium value and unaffected by the deformation 
of the sample. In other words, the equilibrium theory allowed for statistical bond 
fluctuations but not for a dynamical evolution of the parameters characterizing the 
thermodynamic state of the bond network. An obvious starting point for generalizations 
of the model to non-equilibrium situations is therefore to consider the number of closed 
bonds, and therefore also A, as dynamic variables, dependent on the strain- and stress- 
history of the network. We now provide a mean-field description to account for such 
a dynamically evolving bond network. For clarity, we return to the intuitive picture 
underlying the GwLC where the (possibly crosslinker- or molecular-motor-mediated) 
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complex interactions between the polymers are summarized into an effective interaction 
potential for a test segment against the background, as sketched in figure [2] The 
same idea has also been used before in many related situations (e.g. [261 SD 132] )■ 
The generic potential exhibits three essential features: a well-defined bound state, a 
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Figure 2. Sketch of the schematic effective interaction potential for a test segment 
against the background. A polymer segment can either be "bound" by a transient 
crosslinker or a sticky patch on the backbone of another polymer (left potential well) 
or "unbound" and merely confined by the surrounding network to a tube like cage 
(right well). When considering an ensemble of contacts, the fraction of closed bonds 
depends on the free energy difference U between both states. The transition rates fc_ 
and fc+ between the bound and unbound state depend on the barrier heights £ and 
(£ — [/), respectively. An externally applied force can tilt the potential and favour 
binding or unbinding events. 



well-defined unbound state, and a barrier in between (figure [2]). The confining well 
corresponding to the unbound state represents the tube-like caging of the test polymer 
within the surrounding network fi5[ \T7\ Il3] . For ease of notation, we further introduce 
the dimensionless fraction of closed bonds or, bond fraction 

v{t) ^ Ao/A(t) , (3) 

which is simply the ensemble-averaged fraction of sticky contacts that are actually in the 
bound state. The minimum bond distance Aq is typically on the order of Lg, but may be 
somewhat larger in situations where the bonds are mediated by crosslinker molecules or 
by partially sterically inaccessible sticky patches (as e.g. for helical molecules [T2| li^). 

A quantitatively fully consistent way of calculating the dynamics of v would 
involve solving the full Fokker-Planck equation for a Wlc-Wlc contact — a formidable 
program to be pursued elsewhere |15]. Here, for the sake of our qualitative purposes, 
we chose to concentrate onto the physical essence of the discussion, and keep the 
mathematical structure as simple and transparent as possible. We therefore approximate 
the dynamics by a simple exponential relaxation as familiar from the standard example 
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of reacting Brownian particles, conventionally described by Kramers theory with a Bell- 
like force dependence |l6l H?]. Using this simplification and assuming a schematic 
interaction potential as depicted in figure [2} the value of z/ is determined by a competition 
of bond breaking and bond formation with reaction rates fc_ and /c+, respectively. 
Both rates are represented in the usual adiabatic approximation (meaning that the 
equilibration inside the wells is much faster than the barrier crossing and external 
perturbations) by ^B] 

k.To = e-^+(^*-^'')/ , k+To = e-^+^-(^--^*)/, (4) 

where tq is an intrinsic characteristic Brownian time scale for bond breaking and 
formatior][H and / is the force pulling on the bond. Noting that the fraction of open 
bonds is 1 — u, we can then write down the following rate equation for the fraction of 
closed bonds 

The time dependence of u{t) leads via (|3| to an implicit time dependence of the 
GwLC-parameter A{t) and thereby of all observables derived from the GWLC. The 
time- dependent force f{t) in ^ may be thought to result from an externally imposed 
stress protocol or from internal dynamical elements such as molecular motors setting 
the network under dynamic stress. Hence, via an appropriately chosen set of slowly 
changing state parameters f{t), U{t), S{t), . . . the model can in principle accommodate 
for the active biological remodelling in living cells and tissues [SH]. (For a discussion of 
the relation between the microscopic / and the macroscopic stress, see [29].) Note that 
for constant force, the stationary value of u, 

I^stat = {1 + exp[-U + {Xu + Xb)f sta.t])~^ , (6) 

obtained by setting z> = 0, does not depend on £, as it should be (the steady state is 
independent of the transition state). 

3. Results and discussions 

3.1. Vis coelastic properties 

Whenever the bond kinetics can be disregarded (z> = 0), the viscoelastic properties are 
simply those of the bare GwLC [3 123 EO] , which can basically be characterized as 
Wlc short time dynamics followed by a slow, highly stretched logarithmic relaxation 
resembling power-law rheology with a non-universal exponent within the typical 
experimental time windows. Strong perturbations (stress or strain) then give rise to 
a pronounced stiffening due to the steeply nonlinear response of semiflexible polymers 
under elongation [2]. This well established behaviour can be understood as the canvas 

I Strictly speaking the fluctuations in the potential wells at Xb and Xu are characterized by different 
Brownian time scales depending on the width of the potential wells. At the present stage, we do not 
bother to distinguish these time scales nor to fine-tune their numerical values, and identify them, for 
the sake of simplicity, with the entanglement time scale r^^. in our numerical calculations. 
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Figure 3. Dependence of the microrheological "modulus" 5gwlc on the bond fraction 
V for various values of E (brown: E = 0.1, red: £ = 1, blue: £ — 5, green: £ = 10, 
black: £ — oo) at a frequency wo ~ I/IO'T'aq^ and a force / = /ao/2 (taq is the 
relaxation time of a mode of length Aq and /a„ is the Euler buckling force for a rod of 
length Aq); left: normalized real part g'; right: loss angle. 



against which we aim to discern the characteristic mechanical signatures of the bond 
kinetics. 

An observable that will be of particular interest in the following is the complex 
microscopic "modulus" g = 5'gwlo(i^o, Aq/z^, /), introduced in the previous section. It is 



used to determine the linear as well as (via a nonlinear update scheme, see Appendix 



Aj) the nonlinear force response of the system. To understand how the time-dependence 
(|5| of affects this important quantity, it is instructive to first examine the dependence 
of g on u. To this end, we formally consider u temporarily as an independent variable 
instead of determining it from ([s]) . Note that this approach nevertheless makes sense as 
for a fixed set of parameters, u can still take any value between zero and one. This is 
due to the freedom of choosing an initial state, which can be imagined to have evolved 
from the prior deformation history. 

For fixed other parameters, both the real part g' (figure |3^) and the imaginary part 
g" of g increase monotonically with u. We emphasize that g" is not simply proportional 
to g' and therefore the loss angle 6 = aTctan{g" / g') also depends on v (figure |3]d). For 
small z/, the loss angle is large, corresponding to fluid-like behaviour. With increasing 
u, the system becomes more solid-like. Increasing the barrier height S makes the 
dependence oi g on u more pronounced (figure |3]). Conversely, as can be expected, 
the dependence on the barrier height S vanishes with decreasing bond fraction u. Note 
that due to the boundary conditions the limit of a Newtonian fluid [6 = vr/2) is not 
recovered when formally taking the limit of vanishing bond fraction (z/ — )■ 0, A — oo). 
Finally, the potential difference U solely determines the dynamics of u via ([s]), and 
therefore influences the modulus g solely through the equilibrium value for u. 

After these general considerations, we now concentrate on the nonlinear and non- 
equilibrium dynamic response of the extended inelastic GwLC model, which results 
from the coupled relaxation of the viscoelastic polymer network and the transient bond 
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network. 



3.2. Stress-stiffening versus rate- dependent yielding 

A convenient way to characterize the mechanical properties of an inelastic material is a 
force-displacement curve. For a perfectly linear elastic (Hookean) body, it would simply 
consist of a straight line, whereas a perfectly plastic body would feature as a rectangle 
delineated by a yield force fy and an arbitrary plastic strain value. For our qualitative 
purpose, we identify the average transverse displacement of the test polymer segment 



(which is used to determine the force response, see Appendix A) with the reaction 
coordinate x of the schematic potential sketched in figure |2] and we identify the force 
/ entering the reaction rates with the backbone tension of the test polymer. As a 
characteristic length scale for the transverse displacements we use the width {xu — Xt) 
of the effective confinement potential, which is a measure of a typical mean-square 
displacement of the polymer in the unbound state. Using this convention, we now 
consider the effect of a time-symmetric displacement pulse on the evolution of the force 
/(t) and bond fraction i^{t). For technical convenience, we use a Gaussian shape for the 
displacement pulses, but the qualitative conclusions to be drawn are largely independent 
of the precise protocol. The duration of the displacement pulse, which sets the time 
scale for the dynamic response, is used as the unit of time in the following. Here, we 
are not interested in short-time tension-propagation and single-polymer dynamics [H], 
hence a lower bound for pertinent pulse durations is provided by the interaction time 
scale ~ tq. On the other hand, for pulse durations longer than tqc^ the system 
deformation would be quasistatic so that no genuinely dynamic effects of the bond 
kinetics could be observed. 

For a small Gaussian displacement pulse of relative amplitude 0.73 (c.f. figure 



and Appendix A), the shape of the curve predicted by our inelastic GWLC model shows 
all features of a viscoelastic medium (figure |4|3, blue dashed curve). It starts with 
a nearly linear regime for relative displacements < 0.4, where a weak stiffening sets 
in. Due to dissipation in the medium, the path back to the initial state takes its 
course at lower force, which causes a weak hysteresis. This is essentially the viscoelastic 
response that already the bare GwLC model would have predicted. For a larger relative 
deformation amplitude of 2.2, however, the predictions of the bare GwLC and the 
inelastic GwLC diverge. The bare GwLC predicts a strong stiffening (figure |4jD, red 
dotted curve), whereas the full model exhibits an initial stiffening regime followed by a 
pronounced softening (figure |4]d, solid red curve). A very interesting observation is that 
"softening" in this case does not only mean a decrease of the slope of the force-extension 
curve, but that the slope actually becomes negative over an extended parameter region, 
once an operationally defined threshold force fy is exceeded. This "flow state" bears 
more resemblance to plastic flow than to viscoelastic relaxation. The reason for this 
qualitatively new phenomenon is the force-induced bond breaking. Theoretically this 
is best exemplified by the time-dependence of the (experimentally not easily accessible) 
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Figure 4. Transition from a stiffening to a softening response at a rate-dependent, 
operationally defined yield force /y: Gaussian deformation pulses (a) with relative 
displacement amplitudes of 0.73 (dashed) and 2.2 (solid); corresponding force- 
displacement curves and bond fraction evolutions (6, c; dotted curve in b: response of a 
GwLC without bond breaking); order-of-magnitude estimate for the rate-dependence 
of the yield force, ([t]) with const. = 2.2, versus the numerically obtained maximum 
forces of the force-displacement curves for various mean force rates (d); the minimum 
bond fraction (d, inset) depends non-monotonically on the force rate; £ = 10.5, U = 2, 
To = 0.004, {xt ~ Xb)/{xu - Xt) = 0.07, Aq = 22.3. 



fraction z/(t) of closed bonds (figure |4p). A glance at the bond fraction during the weak 
deformation scenario (figure |4|:;, blue dashed curve) reveals how the limit of a bare GwLC 
is recovered, namely whenever the deformation is not sufficiently violent and persistent 
to significantly decrease the fraction of closed bonds. For the large deformation scenario, 
in contrast, the bond fraction stays only initially constant (figure |4|3, solid red curve). As 
a consequence of the large strain the force rises steeply, as can be seen from the strong 
stiffening in figure |4]d (sohd red curve). At the yield force, the bond fraction suddenly 
decreases to nearly half of its initial value during a very short time. The decrease in the 
bond fraction is accompanied by a somewhat slower drop in the force, which is reflected 
in the sudden softening of the force-displacement curve. The bond fraction eventually 
recovers on a much slower time scale, which is roughly given by tq exp{S — U). Note that 
even though the effects of the inelastic response dominate the stress-strain curve, the 
viscoelastic relaxations from the underlying GwLC model are still present. They could 
hardly be disentangled from the inelastic contributions, though, without an underlying 
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faithful model of the viscoelastic response at hand. 

In summary, we observe a competition between force-amplitude dependent stress- 
stiffening and force-rate dependent yielding events. If the backbone force stays much 
smaller than the yield force, the adaptation of the bond network requires an adiabatically 
long time, on the order of the bond lifetime tqc^ , so that plain viscoelastic GWLC 
behaviour is observed on the pulse time scale. In contrast, if the backbone tension 
/ reaches the yield force /y, the bond fraction declines sharply, thereby switching the 
response from the Wlc / GwLC-typical stiffening to a pronounced softening. 

The rate dependence of the yield force fy can be estimated by setting the time it 
takes to reach the yield force, fy/ f, equal to the force-dependent time scale of bond 
opening, kZ^{f), from 

fyT^ ~ kz\fy) => {xt - Xb)fy ^ LW (^const. x {xt - Xb)Toe^ f^ ■ (7) 

Here, LW(x) is the positive real branch of the Lambert W function. In figure |4]i, this 
estimate is compared with results from numerical evaluations for Gaussian displacement 
pulses at different average rates. Apart from the numerical errors, the slight deviations 
from the estimate may be attributed to the fact that the force rate is not constant for 
the Gaussian protocol. They can be mostly eliminated by using force ramp protocols of 
various slopes instead of the Gaussian displacement pulses. For not too low rates, the 
rate dependence can be approximated by the even simpler relation 

{xt - Xb)fy ^{S + ln f'To + const.) , (8) 

where the force rate has to be normalized by a suitable force scale. 

The minimum bond fraction reached during the application of the pulse, which 
we interpret as a measure of the degree of fiuidization, depends non-monotonically on 
the rate (figure |4]i, inset). Qualitatively, this can be understood by noticing that two 
different factors influence the fiuidization, namely the maximum force attained during 
deformation and the time over which the force is applied. The maximum force is simply 
the yield force /y, with the rate dependence in ([T]) and figure |4]i, while the reciprocal 
rate itself sets the time scale. For high rates {xu — Xt)f > 1 (in our example) the 
rate-dependence of the maximum force wins and the minimum v reached decreases with 
increasing rate. For low rates, where the rate-dependence of the maximum force is 
much weaker (figure |4]i) the minimum u decreases with increasing pulse duration, viz. 
decreasing rate. Note that for slow pulses (low rates), the bond fraction after the pulse 
may be quite different from the minimum z/, as significant recovery may already take 
place during the pulse. 

3.3. Transient remodelling and cyclic softening 

The substantial changes in the material properties accompanying bond breaking can 
be exemplified by monitoring the linear elastic modulus g'{u}o) (measured at a fixed 
frequency Uq) in response to a strain pulse (figure |5^). Apart from the usual 
WLc/GwLC-typical stress-stiffening below the threshold force fy, the modulus is seen 
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to drop below the value it had before the deformation pulse, where it apparently 
saturates (remember that the deformation vanishes for t > 1 and that the recovery 
takes roughly a time tqc^"^ ^ 1). We thus observe what we call a "passive", "physical" 
remodelling of the bond network as opposed to the "active" , "biological" remodelling of 
the cytoskeleton of living cells in response to external stimuli. These passive remodelling 
effects have recently been observed for human airway smooth muscle cells [31] . (A more 
quantitative discussion of the relation to the experimental observations will be given 
elsewhere.) The fact that the deformation pulse leads to a decrease in g' and an increase 
in the loss angle S (figure [sj^) suggests the notion of fluidization [31] . 



a b 




t t 



Figure 5. The internal force / and tlie bond fraction v largely determine the storage 
modulus g' and the loss angle 5: normalized storage modulus g'^ (a) and loss angle 5 (c) 
during displacement pulses of large and small amplitudes, 2.2 (solid), 0.73 (dashed); 
the stronger deformation pulse causes a strong and persistent overshoot of the loss 
angle, indicating fluidization] slow recovery after the pulse, (6, d); parameters as in 
figure |4) i:jo = 3. 



The passive remodelling described here is a transient phenomenon, because, after 
cessation of the external perturbations, the bond fraction will eventually recover its 
equilibrium value. This indicates that also the change in the system properties is a 
transient effect, which is demonstrated by the recovery from fluidization in figure [5] b 
& d. Thus, while the fluidization resembles a plastic process on short time scales, on 
long times scales, the phenomenology is more similar to pseudo- and superelasticity as 
observed in shape-memory alloys |49j . 
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To demonstrate that the transient passive remodelhng also affects the nonhnear 
material properties, we consider a series of three pulses, next (figure [6]). The force 
response to such a protocol is depicted in figure [6]d. The force-displacement curve 
for the second stretching (second left branch of the solid curve) is less steep than for 
the first stretching (very left branch of the solid curve) and the yield force is lower. 
This indicates a cyclic softening, or viscoelastic shake-down effect, closely related to 
the fiuidization of the network by strain. The strength of the shake-down depends 
on the fraction of transiently broken bonds, and hence on the rate as well as on the 
amplitude of the applied deformation (c.f. figure [6]:, solid red curve). Upon repeated 
application of deformation pulses the force-displacement curves settle on a limit-curve 
that is essentially preconditioned by the initial deformation pulse and the inelastic work 
it performed on the sample. For more gentle protocols that only break a smaller fraction 
of bonds, the initial fiuidization would be not as pronounced and one would obtain a 
gradual shake-down, which converges to a limit curve only after many cycles. 




Figure 6. Protocol (a) and response {b, c) for three subsequent Gaussian strain pulses; 
force-displacement curves (6) for a single minimum interaction distance (Aq — 22.3, 
solid) and for a distribution of Aq according to [40; with mean Aq = 22.3 (dashed); the 
bond fraction i' (c) for Aq = 0.5 Aq (dash-dotted), Aq (solid), and 2Ao (dashed); other 
parameters as in figure |4] 
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3.4- Introducing multiple length scales 

So far, we have assumed that there is one well-defined characteristic length scale Aq 
for the polymer interactions, which is on the order of the entanglement length Lg and 
interpreted as the minimum average contour distance between adjacent bonds of the test 
polymer with the background network. This is clearly a mean-field assumption. Recent 
combined experimental and theoretical studies have established that the local tube 
diameter and entanglement length in pure semifiexible polymer solutions actually exhibit 
a skewed leptocurtic distribution with broad tails [IQ]. To include this information 
in the above analysis is not straightforward, since the elastic interactions between 
regions of different entanglement length are not known, a priori. For the sake of a first 
qualitative estimate, the simplest procedure seems to be to average the above results over 
a distribution of Aq, corresponding to a parallel connection of independent entanglement 
elements. Qualitatively, this renders the abrupt transition from stiffening to softening 
somewhat smoother (c.f. figure [6]d, dashed curve). Moreover, the initial stiffening is 
slightly more pronounced. This is due to the contributions of Aq < Aq, which exhibit 
both stronger stiffening and fluidization. Nevertheless, all qualitative features — like 
stiffening and remodelling effects and a "flow state" indicating fluidization or inelastic 
shakedown — can still be well discerned. This is consistent with the assumption that 
pure semiflexible polymer solutions are at least qualitatively well described by a single 
entanglement length [T71 [T8| [22| [271 [50] . 

4. Conclusions and outlook 

We have presented a theoretical framework for a polymer-based description of the 
inelastic nonlinear mechanical behaviour of sticky biopolymer networks. We represented 
the polymer network on a mean-fleld level by a test polymer described by the 
phenomenologically highly successful GWLC model. Beyond the equilibrium statistical 
fluctuations captured by the original model, we additionally allowed for a dynamically 
evolving thermodynamic state variable A{t) characterizing the network of thermo- 
reversible sticky bonds, which is updated dynamically upon bond breaking according 
to an appropriate rate equation. Within our approximate treatment we could obtain 
a number of robust and qualitatively interesting results, which are not sensitive to the 
precise parameter values chosen, nor to the choice of the transverse rather than the 
longitudinal susceptibility in ([T]). For sufficiently strong deformations, changes in the 
mean fraction of closed bonds, refiected in A(t), can infiuence material properties to a 
degree at which the material behaves qualitatively different from the usual viscoelastic 
paradigm. In particular, our theory predicts a pronounced fluidization response, which 
develops upon strong deformations on top of the intrinsic viscoelastic stiffening response 
provided by the individual polymers. After cessation of loading, the system slowly 
recovers its initial state. These observations are in qualitative agreement with recently 
published experimental results obtained for living cells [34J, and a more quantitative 
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comparison with dedicated measurement results for cells and in vitro biopolymer 
networks will be the subject of future work. We found the yield force for the onset 
of the fiuidization to be sensitive to the deformation rate. Moreover, the fluidization 
response was shown to be accompanied by a cyclic softening or shake-down effect. 
Taking into account the spatial heterogeneity of biopolymer solutions by a distribution 
of entanglement lengths leads to a smoothing of the force-displacement curves without 
affecting their qualitative characteristics. One may still expect qualitatively new effects 
in situations with unusually broad distributions of entanglement lengths, such as for 
strongly heterogeneous (e.g. phase separating) networks. 

A parameter that was found to be very important for cells [51] but has not 
been discussed much in this contribution is the prestressing force fo (see appendix 



Appendix A), which is present in adhering cells even in the absence of external driving. 



Experimentally, it was established that increasing prestress is correlated with a higher 
stiffness of adhering cells [51] which in turn is correlated with a higher stiffness of 
the substrate [52]. When naively treating the prestressing force just as an additive 
contribution to the overall force, our model predicts the opposite: the additional force 
breaks more bonds and the network becomes softer and more fluid-like. To reconcile 
these apparently contradicting trends, one can appeal to the notion of homeostasis, 
which essentially amounts to postulating that the cell actively adapts its structure such 
that a certain set of thermodynamic state variables remain in a physiologically sensible 
range [39] . This basically implies that the cell will respond to stiff substrates or persistent 
external stresses by a biological remodelling that corresponds to an increase of S and 
U, and possibly to a decrease of Aq. By adapting also the internal prestress fo, the 
cell may then avoid, apart from the structural collapse, an equally undesirable loss of 
flexibility. We note, however, that internal stresses are indeed observed to disrupt the 
cytoskeleton, as implied by the simple physical theory presented here, if they are not 
permanently balanced by a substrate [53] . 

Even though the nonlinear effects presented in this contribution depend crucially on 
the bond kinetics, it should not be overlooked that the GwLC is an equally indispensable 
ingredient of the complete theory. First, it provides the constitutive relation that gives 
the force in response to an infinitesimal displacement, which in turn governs the bond 
kinetics. Secondly, the time-dependent linear susceptibility strongly filters the dynamic 
remodelling of the bond network in practical applications, in which one rarely has access 
to the microstate of the underlying bond network. In summary, the extended model 
thus integrates experimentally confirmed features of the GwLC response such as slow 
relaxation, power-law rheology, viscoelastic hysteresis, and shear stiffening with a simple 
bond kinetics scheme, which allows less intuitive complex nonlinear physical remodelling 
effects like fiuidization and inelastic shake-down to be addressed. 

Finally, we want to point out that it is straightforward to extend the above analysis 
in a natural way to account for irreversible plastic deformations [15]. To this end, one 
has to allow for the possibility that the transiently broken bonds reform somewhere else 
than at their original sites (as always assumed in the foregoing discussion) after strong 
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deformations with finite residual displacement. This can be included by accommodating 
an additional term acting; ClS db "slip rate" in ^ [15]. 
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Appendix A. Obtaining force-displacement curves and solving for u 

The deformation pulse protocol used in this paper is 



it 



r 



a 



r(t)=7oexp(-^-^). (A.i; 



To obtain the full nonlinear response of the (extended) GwLC to the displacement 



protocol (A.I), we make use of the superposition principle: We know that the Fourier 
transform of QQ^i^ci!^ ^ f) / i}^) is the linear response to a small displacement step. 
After decomposing a finite displacement into infinitesimal displacement steps and a 
partial integration, we write 

/oo POO J 

dt' / -=^?GwLc(u;,A(t),/GwLc(t))e-^"(*-*')r(t'). (A.2) 
-oo J — oo V 



Note that the right hand side of (A.2) depends on /gwlc(^), rendering it a highly 



nonlinear implicit equation. Using (A.I) and integrating by parts, we obtain 



/GwLc(t) = v^7o^ / c?a;e-*"(*-")(7GwLc(u;,A(t),/GwLc(t))e-^. (A.3) 







The next step is to identify a connection between the force in ([s]) and /GwLc(t). A 
finite prestressing force /o is introduced mainly for technical reasons, namely to avoid 
the unphysical region of negative (i.e. compressive) backbone stress, which would buckle 
the polymers. While the physics of the prestressed network can essentially be mapped 
back to that of a network without prestress by a renormalization of the parameters S 
and U, prestress may also be seen as an important feature: indeed, the cytoskeleton of 
adhered cells is well known to be under permanent prestress, and suspended cells seem 
prone to spontaneous shape oscillations that could be indicative of a propensity of the 
cells to set themselves under prestress [3S] . In the context and on the longer time scales 
of the biological remodelling of the cytoskeleton it would therefore make sense to think 
of /o as a dynamic force generated by molecular motors and polymerization forces. For 
the following, however, we take the prestress to be constant so that the force entering 
(§ is 

/(t) = /o + /GwLc(t). (A.4) 
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We now face the problem that (A. 2) is an imphcit equation for /gwlc(^), which depends 
on z/(t), which in turn depends via f{t) on /gwlc(^)- We therefore use a two-step Euler 
scheme: As initial values for / and we choose the prestressing force, f{t = 0) = /o, 
and the steady-state value under the prestressing force, z/(t = 0) = z/(t — )■ oo|/ = /q), 
respectively. We then choose a sufficiently small time step At and apply the following 
iterative rule: 

f{kM) ^ /o + fawLc (kAtlfiik - l)At, u{{k - l)At)) (A.5) 
i>{kAt) ^ i>{kAt\f{kAt)), (A.6) 
where k G N+. In the limit At — t- 0, this procedure converges to the exact solution. 
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